# -*- coding: utf-8 -*-
"""
PBAR_Zspec_DetermineBestAlgorithmAllCounts2ndSet.py

Created on Wed Apr 16 09:22:57 2014

@author: jkwong
"""

#import csv
#import matplotlib.pyplot as plt
import os, cPickle, copy
import numpy as np
import matplotlib.pyplot as plt
#import numpy.matlib
#import datetime
#import time
import PBAR_Zspec
reload(PBAR_Zspec)
from sklearn.feature_selection import SelectKBest, chi2, f_classif, f_oneway
from sklearn.linear_model import LogisticRegression
from sklearn.lda import LDA
from sklearn import cross_validation
from sklearn.metrics import precision_recall_fscore_support
from sklearn.preprocessing import scale
reload(PBAR_Zspec)

def gauss_function(x, a, x0, sigma):
    return a*np.exp(-(x-x0)**2/(2*sigma**2))

# Set useful plot variables
plotColors = ['r', 'b', 'g', 'm', 'c', 'y', 'k'] * 10
lineStyles = ['-', '-.', ':', '_', '|'] *  10
markerTypes = ['.', 'o', 'v', '^', '<', '>', '1', '2', '3', '4', 's', 'p', '*', 'h', 'H', '+', 'x', 'D', 'd']
markerTypes = markerTypes * 2


setNum = 3 # 2a

# First set
if setNum == 1:
    basepath = r'C:\Users\jkwong\Documents\Work\PBAR\data'
    infoFilename = os.path.join(basepath,'datasetSummaryOLD.txt')
# Second set
elif setNum == 2 or setNum == 3 or setNum == 4:
    # SECOND SET
    basepath = r'C:\Users\jkwong\Documents\Work\PBAR\data4\Mar-files'
    infoFilename = os.path.join(basepath,'DatasetSummary2ndSet.txt')

fullFilename = os.path.join(basepath, 'zspec%dSet.dat' %setNum)
with open(fullFilename ,'tb') as fid:
    print('Writing %s' %fullFilename)
    datAll = cPickle.load(fid)




# try one of the classifers on the whole set
countRangeIndex = 16
numberFeaturesIndex = 4
#trialIndex = 0
lda = classifiersWholeAllAll[countRangeIndex][numberFeaturesIndex]['lda']
topFeaturesIndices = topFeaturesIndicesAllAll[countRangeIndex][numberFeaturesIndex]
countRange = countRangesList[countRangeIndex]
cutCount = (countMatrix > countRange[0]) & (countMatrix < countRange[1])

# define feature target matrices
X = statsMatrix[cutCount,:][:,:-1]
X_mean = X.mean(0)
X_std = X.std(0)
#X_scaled = scale(X)
X_scaled = copy.copy(X)
y = statsMatrix[cutCount,:][:,-1]

#X_All_scaled =( statsMatrix[:,:-1] - np.tile(X_mean, (statsMatrix.shape[0], 1)) )/ np.tile(X_std, (statsMatrix.shape[0], 1))

X_All_scaled = copy.copy(statsMatrix[:,:-1])

y_All = statsMatrix[:,-1]
X_All_transformed = lda.transform(X_All_scaled[:,topFeaturesIndices])
X_decision_function = lda._decision_function(X_All_scaled[:,topFeaturesIndices])
pred= lda.predict(X_All_scaled[:,topFeaturesIndices])


X_All_transformed2 = np.dot(X_All_scaled[:,topFeaturesIndices] - lda.xbar_, lda.scalings_)

n_comp = X_All_transformed2.shape[1] if lda.n_components is None else lda.n_components
X_All_transformed3 = np.dot(X_All_transformed2, lda.coef_[:n_comp].T)
X_All_transformed4 = np.dot(X_All_transformed2, lda.coef_.T) + lda.intercept_
X_decision_function = lda._decision_function(X_All_scaled[:,topFeaturesIndices])


plt.figure()
plt.grid()
cut = y_All == 0
plt.plot(countMatrix[cut], X_decision_function[cut,0] - X_decision_function[cut,1], '.b')
cut = y_All == 1
plt.plot(countMatrix[cut], X_decision_function[cut,0] - X_decision_function[cut,1], '.r')
#cut = pred == y_All
#plt.plot(countMatrix[cut], X_decision_function[cut,0] - X_decision_function[cut,1], 'om', alpha = 0.1, markersize = 13)

cut = pred != y_All
plt.plot(countMatrix[cut], X_decision_function[cut,0] - X_decision_function[cut,1], 'sc', alpha = 0.1, markersize = 13)
title('Trained on %d, %d' %(countRange[0], countRange[1]) )



plt.figure()
plt.grid()
cut = y_All == 0
plt.plot(countMatrix[cut], X_All_transformed2[cut], '.b')
cut = y_All == 1
plt.plot(countMatrix[cut], X_All_transformed2[cut], '.r')

cut = pred == 1
plt.plot(countMatrix[cut], X_All_transformed2[cut], 'om')


plt.ylabel('Discrim')
plt.xlabel('Counts')


plt.figure()
plt.grid()
cut = y_All == 0
plt.plot(countMatrix[cut], X_decision_function[cut,0], '.b')
cut = y_All == 1
plt.plot(countMatrix[cut], X_decision_function[cut,0], '.r')
cut = pred == 1
plt.plot(countMatrix[cut], X_decision_function[cut,0], 'om', alpha = 0.1, markersize = 10)
